i = 0
rf.results <- list()
for (year in start.year:end.year) {
  i = i + 1
  print(year)
  # Generate dependent variable
  Violence.mean <-  aggregate(x = as.matrix(dta[dta[,t.var]< year,
                                                dv]),
                              by = list(dta[dta[,t.var]< year,
                                            loc.var]),
                              FUN = mean)
  colnames(Violence.mean) <- c(loc.var,
                               "dv_mean")
  Violence.new <- merge(as.matrix(dta[dta[,t.var]< year,
                                      c(loc.var,t.var,dv)]),
                        Violence.mean,
                        by = loc.var)
  Violence.past <- Violence.new[,dv]-Violence.new$dv_mean
  # Fit Random Forest
  set.seed(i*123)
  rf.results[[i]] <- randomForest(x = as.matrix(dta[dta[,t.var]< year,
                                                         rhs]),
                           y = as.matrix(Violence.past),
                           ntree = rf_demean_params$ntree,
                           nodesize = rf_demean_params$nodesize,
                           mtry=rf_demean_params$mtry,
                           type = rf_demean_params$type,
                           do.trace = rf_demean_params$do.trace,
                           importance = rf_demean_params$importance
  )
  rf.results[[i]]$optim.size <- 15
  rf.results[[i]]$fit.oos <- predict(rf.results[[i]],
                                        newdata = as.matrix(dta[dta[,t.var] == year,
                                                                  rhs]),
                                     type = "response")
  Violence.new <- merge(as.matrix(dta[dta[,t.var]== year,
                                           c(loc.var,t.var,dv)]),
                        Violence.mean,
                        by = loc.var)
  rf.results[[i]]$actual.oos <- Violence.new[,dv]-Violence.new$dv_mean
  # Display time
  hours <- floor((proc.time()[3]-start.time)/3600)
  mins <- floor((proc.time()[3]-start.time)/60) - hours*60
  secs <- floor((proc.time()[3]-start.time)) - hours*3600 - mins*60
  
  print(paste(hours,
              "h",
              mins,
              "m",
              secs,
              "s elapsed",
              sep = ""))
}
save(rf.results,file=paste(modeldir,"/",
                              country,
                              "_rf_",
                              "demean",
                              "_",
                              rhs.group,fileext,
                              ".RData",
                              sep = ""))